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ABSTRACT 

Let  {X(t )  :  t  >  0}  be  a  real-valued  etochastlc  process  and  set 
a  «  j"  X(t)  G(dt),  where  G  is  a  (non-random)  distribution  function.  If 
support  of  G  is  large,  standard  Monte  Carlo  techniques  for  estimating  a 
are  inefficient,  since  X  must  be  simulated  over  the  entire  support  of  G 
To  avoid  this  difficulty,  randomization  schemes  are  derived  that  require 
simulation  of  X  over  random  subsets  of  the  support  of  G.  Large-sample 
behavior  of  randomized  estimators  is  studied  in  detail.  Some  variance 
reduction  techniques  are  also  presented. 


AMS  (MOS)  Subject  Classifications:  60G57,  68J05 

Key  Words:  Simulation,  variance  reduction,  random  measures 

Work  Unit  Number  5  -  Optimization  and  Large  Scale  Systems 


Sponsored  by  the  Uh)ted  States  Army  under  Contract  No.  DAAG29-80-C-0041 


SIGNIFICANCE  AND  EXPLANATION 


-j  V"'  C  *■  «■-  f-  y  L 

Consider*^  stochastic  system  for  which  one  needs  an  estimate  of  the 
expected  discounted  cost  over  the  infinite  horizon.  Standard  Monte  Carlo 
procedures  do  not  apply  since  the  parameter  to  be  estimated  involves  values  of 
the  process  over  an  infinitely  long  time  interval/  In  T-his  paper,  we  presents 
Monte  Carlo  estimation  techniques,  based  on  randomization,  that  can  be  used  in 
the  above  setting.  The  techniques  developed  turn  out  to  be  more  efficient 
than  the  standard  approach,  even  when  the  parameter  to  be  estimated  cumulates 
costs  over  a  finite  time  interval. 


\ 
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The  responsibility  for  the  wording  and  views  expressed  in  this  descriptive 
summary  lies  with  MRC ,  and  not  with  the  author  of  this  report. 


RANDOMIZED  ESTIMATORS  FOR  TIME  INTEGRALS 
Peter  w.  Glynn 

1.  Introduction 

Let  {x(t):t>0}  be  a  real-valued  stochastic  process  representing  the 
output  of  a  simulation,  and  consider  a  time  integral  of  the  form 

(1.1)  1  -  J"  X(t)  G(dt)  , 

where  G  is  a  (deterministic)  non-decreasing  function.  Our  objective,  in 
this  paper,  is  to  study  Monte  Carlo  techniques  for  estimating  the  parameter 
a  -  El. 

The  time  integral  (1.1)  includes  several  commonly  studied  performance 
criteria. 

(1.2) .  Example.  If  G(t)  «  1  -  e  (o>0),  then  I  corresponds  to 
discounting  X(t)  at  rate  a  over  an  infinite  horizon.  Such  time  integrals 
occur  frequently  in  inventory  models. 

(1.3) .  Example.  If  G(t)  ■  min(t/T,  1),  then  I  is  the  average  of  X(t) 
over  the  interval  [0,T].  Such  averages  are  often  of  interest,  in  a  queueing 
context. 

(1.4) .  Example.  Let  {Y(t)  i  t  >  0)  be  a  stochastic  process  and  f(t,y)  be 
a  real-valued  performance  measure  which  may  depend  explicitly  on  t.  By 
setting  X(t)  ■  f(t,Y(t)),  criteria  of  the  form 

I  ■  Jq  f ( t ,Y(t) )  G(dt) 

can  be  incorporated  as  a  special  case  of  (1.1). 

(1.5) .  Example.  Grassman  (1982)  has  recently  developed  a  Monte  Carlo 
technique  for  estimating  Ef(Z(t)),  where  Z(t)  is  a  uniformizable  Markov 
jump  process.  The  idea  is  to  represent  Z(t)  as  Y(N(t)),  where  N(t)  is  a 


Poisson  process  with  rate  X  (say),  and  {Y(k)  :  k  >  0}  is  an  independent 


discrete-tine  Markov  chain.  Thus, 

(1.6) 


Ef (Z(t) )  -  l  Ef(Y(k))e"Xt(Xt)k/kl  . 
k-0 


Setting  X(t)  -  f(Y([t]))  ( [t]  -  greatest  integer  less  than  or  equal  to  t), 
and  G(dt)  -  Poisson  measure,  we  see  that  the  representation  (1.6)  is  a 
special  case  of  (1.1). 

In  Section  2,  we  shall  briefly  discuss  the  direct  method  for  estimating 
a;  this  involves  generation  of  independent  variates  1^,  I^,...,  each  having 
the  distribution  of  X.  The  parameter  a  is  then  estimated  by  I(n),  where 


n 


1(n)  *  n  E  \  * 

n  k-1  K 


The  difficulty  with  the  direct  method  is  that  if  T  ■  sup{t  :  G(t)  <  G(<*)} 


(G(*»)  “  limG(t))  is  large,  then  generation  of  variates  is  expensive.  Thus, 

t-M* 


in  Section  3,  a  general  framework  for  randomized  estimation  of  a  *  El  is 
presented. 

(1.7).  Definition.  R(n)  is  a  randomized  I(n)  -  estimator  if  there  exists  a 
o-field  G  such  that 


e{r( n)  |  G)  *  I(n)  . 

The  definition  of  conditional  expectation  implies  that  ER( n)  =  EI(n)  * 
a,  justifying  the  description  of  R(n)  as  an  I(n)  -  estimator.  Some 
authors  refer  to  methods  based  on  randomized  estimation  as  conditional  Monte 
Carlo  procedures  (see,  for  example,  Rubinstein  (1981),  p.  141).  However,  we 
prefer  to  reserve  the  term  "conditional  Monte  Carlo"  for  the  "converse"  to  a 
randomized  estimator. 

(1.8).  Definition.  I(n)  is  a  conditional  Monte  Carlo  R(n)  -  estimator  if 
and  only  if  R(n)  is  a  randomized  I(n)  -  estimator. 


The  definition  (1.8)  of  the  term  "conditional  Monte  Carlo"  is  consistent 
with  its  usage  in  several  recent  books;  see,  for  example,  Bratley,  Fox,  and 
Schrage  (1983)  or  Law  and  Kelton  (1982).  The  following  proposition  is  a  well- 
known  property  of  conditional  expectation  (Burrill  (1972),  p.  392). 

(1.9).  Proposition.  Suppose  that  E|R(n)|  <  •  and  that  R(n)  is  a 
randomized  I(n)  -  estimator.  Then, 

var(R(n))  >  var(I(n)) 

Proposition  1.9  states  that  a  randomized  I(n)  -  estimator  has  larger  mean 
square  error  (MSB)  than  I(n).  For  this  reason,  the  conditional  Monte  Carlo 
estimator  I(n)  is  to  be  preferred  in  the  case  that  generating  ft(n) 
requires  the  same  effort  as  simulating  I(n).  However,  it  turns  out  that  in 
our  time  integral  setting,  the  time  required  to  generate  a  randomized  X(n)  - 
estimator  will  often  be  smaller  than  that  required  to  simulate  T( n) ;  this 
property  can  offset  the  MSB  advantage  of  I(n).  This  is  the  theme  of  Sections 
3  through  6.  After  developing  a  general  framework  for  randomized  estimation 
in  Section  3,  three  specific  randomized  estimation  algorithms  are  studied  in 
Section  4  through  6. 

2.  The  Direct  Method 

Development  of  estimation  theory  for  the  parameter  a  requires  some 
assumptions  on  the  simulation.  Let  (ft,F,P)  be  the  probability  space  which 
supports  our  simulation.  We  assume  that: 

A1 .  there  exists  a  family  of  processes  {(X^S^)  :  k  >  1}  such  that 
^  :  R+  x  n  ♦  R  and  :  *+  x  0  ♦  *+,  where  R+  =  [0,«) 

A2.  the  processes  X^(*,w)  possess  left  limits  and  are  right  continuous. 


for  each  u>  E  0 


(2.2) .  Theorem .  Assume  A1  -  A6.  Then,  I(t)  ♦  a  a.s.  as  t  +  •  if 

P<Sk(T)  <  ->  -  1. 

Proof »  The  strong  law  of  large  numbers  guarantees  that  I(n)  ♦  a  a.s.  By  A3 

and  A5,  P{Sk(T)  >  0}  >  0  so  N(t)  <  ®  a.s.  for  all  t.  Furthermore,  the 

assumption  P{Sk(T)  -  ■  0  assures  that  N(t)  ♦  »  a.s.  (see  Cinlar 

(1975),  p.  290),  from  which  the  result  follows.  || 

It  is  worth  observing  that  if  P{sk(T)  -  •}  >  0,  then  I(t)  does  not, 

in  general,  converge  to  o  (see  Example  1.2).  we  will  also  be  interested  in 

rates  of  convergence  for  our  estimators.  For  such  results,  we  require  a 

further  moment  assumption. 

A7.  o2  -  E(I  -  a)2  <  «  . 

k 

(2.3) .  Theorem.  Assume  A1  -  A7.  if  p{sk(T)  <  «•}  ■  1,  then 

(2.4)  15  a(t)  •  1 1 ( t )  -  a|  -  a(ES  (T))1^  a.s. 

* 

1/> 

where  a(t)  A  (t/2  log  log  t)  2  and  0.®  A  0. 

Proof.  The  Hartman-Wlntner  form  of  the  law  of  the  iterated  logarithm  implies 
that 

(2.5)  lim  a(n)|l(n)  -  a|  =  a  a.s. 

n-»® 

Since  N(t)  ♦  ®  a.s.,  (2.5)  yields 

(2.6)  lim  a(N(t))|l(t)  -  a|  ®  a  a.s. 

t+® 

But  N(t)/t  *  1/ESk(T)  a. 8.  (see  (6],  p.  290),  which  in  turn  implies  that 

(2.7)  a(N(t))/a(t)  ♦  1/(ESk (T) ft  a.s. 

Relation  (2.4)  follows  immediately  from  (2.6)  and  (2.7).  |  | 

Confidence  intervals  for  a  can  also  be  constructed  from  the  direct 
estimator  I(t).  The  key  tool  is  a  central  limit  theorem  (CLT),  which  is 
valid  under  slightly  stronger  assumptions  than  Theorem  2.3. 

(2.8) .  Theorem.  Assume  Al  -  A7.  if  ES.  (T)  <  ®,  then 


(2.9)  t^  (I(t)  -  a)  — >  O(ESk(T))/2N(0,1  ) 

»•  t  ♦  where  N(0,1)  is  a  mean  zero  unit  variance  r.v.  and  ==>  denotes 
weak  convergence. 

Proof.  From  the  CLT  for  i.i.d.r.v.'s, 

1a  _ 

n'2  (I (n)  -  a)  — >  o  N ( 0, 1 ) 

as  n  ♦  •.  Let  t^  be  an  arbitrary  sequence  converging  to  infinity.  Since 
mv/tk  *  a.s.,  one  can  apply  Theorem  7.3.2  of  Chung  (1974)  to 

conclude  that 

1a 

N(tkr2  (I(tk)  -  a)  — >  0  N (0 , 1  ) 

as  k  ♦  •.  The  converging-together  lemma  (Billingsley  (1968),  p.  25)  then 
yields 

1/_  U 

tk2  (I(tk)  -  a)  ==>  o(ESk(T))/2N(0,1) 

as  k  ♦  Since  {tkJ  was  arbitrary,  we  obtain  (2.9)  ([2],  p.  16).  || 

If  Zg  solves  p{n (0,1)  <  z^}  ■  1  -  6/2,  Theorem  2.8  proves  that  the 
random  interval 

(I(t)  -  zfiv/2  (t)/t/2  ,  I(t)  +  z^v^tj/t^ 

is  an  approximate  100(1-6)%  confidence  for  o,  where  v(t)  is  a  con- 

2 

sistent  estimator  for  o  ES,  (T). 

k 

Before  concluding  this  section,  we  consider  a  special  case  of  Example 
1.4.  Suppose  that  the  process  Y(t)  considered  there  is  a  stationary  process 
on  state  space  S.  Then, 

(2.10)  a  -  El  =  J“  Ef (t,Y(t) )  G(dt) 

=  j"  js  f(t,y)P{Y(t)  e  dy}  G(dt) 

=  Js  Jq  f(t,y)  G(dt)ir (dy) 


where  *(•)  =  P{Y(0)  e  •}.  In  general,  if  the  process  Y(t)  can  be 


simulated,  then  the  distribution  if  is  known  explicitly  so  that  a  can  be 
calculated  analytically  from  (2.10).  Thus,  the  interest  in  Monte  Carlo 
estimation  of  a  occurs  when  Y(»)  is  non-stationary.  In  the  simulation 
literature,  such  non-stationary  estimation  problems  are  referred  to  as 
transient  simulations. 


3.  A  General  Framework  for  Randomized  Estimation  of  Time  Integrals . 

In  this  section,  we  assume  that: 

A8.  there  exists  a  sequence  of  processes  {H^  :  k  >  1}  such  that 
:  R  x  ft  ♦  R+  and  which  satisfy: 

i.)  (•,<!>)  is  non-decreasing  and  right  continuous,  for  each  u  e  (2 

ii.)  {Hj'  :  k  >  1}  is  a  sequence  of  i.i.d.  random  elements,  which  is 
independent  of  the  collection  {(X^,Sk)  :  k  >  1} 

iii.)  EH^t)  -  G(t)  for  all  t  e  R. 

(3.1).  Lemma.  Under  A5  and  A8,  8^(0-)  *  0  a.s.  and  EH^(°*)  »  1. 

Proof.  For  t  <  0,  Hk(t)  is  non-negative  and  satisfies  EHk(t)  =  G^t)  «*  0, 
so  that  Hfc(t)  *  0  a.s.,  proving  that  Hj^O-)  *  0  a.s.  For  the  finiteness  of 
EH^  («• ) ,  observe  that 
convergence, 

(<*)  *  lim  EH.  (n)  ■  lim  G(n)  *  1  I  I 

n+«  n>» 

Our  randomized  I(n)  -  estimator  will  be  based  on 


H^(n)  increases  to  H^(®) .  Hence,  by  monotone 


(3.2) 


V“> 


f  rt  » t  \ 

\  0  ;  (d  6  Aj^ 


where  Ak 


tJ0Uk(t)|Hk(dt)  *  •}.  Observe  that  for  each  u>,  R^(<»»)  is  well- 


defined  by  formula  (3.2),  since  X^  is  product  measurable,  and  satisfies 

A8  i.).  Several  later  arguments  will  require  the  following  approximation 


result.  Let 


X^ft)  ■  mintx^ftjjm} 

X*  (t)  -  min{X  (2”n(  j+1 ) )  ,m} 

KB  K 

for  j  <  2nt  <  j+1. 

(3.3) .  Proposition.  Assume  A1  -  A6,  A8,  and  suppose  that  X^  is  non¬ 
negative  .  Then 

(3.4)  lim  lim  j“  xj^(t)Hk(dt)  *  j"  Xk(t)H  (dt)  a.s. 

m-M*  n-^» 

Proof .  First,  it  is  easily  seen  that  is  right  continuous.  The 

right  continuity  of  X^  implies  that 


lim  X“  (t) 
km 

n-H» 


Xkm(t> 


for  all  t  >  0.  Now,  Hk{<»)  <  ®°  a*s*»  so  that  we  can  apply  a.s.  the 
bounded  convergence  theorem  to  conclude  that 

lim  r0  Km(t)Hk(dt)  =  Jo  Xkm(t)Hk(dt>  a*8* 


Now,  apply  the  monotone  convergence  theorem  to  complete  the  proof .  | | 

Since  the  integrals  on  the  left-hand  side  of  (3.4)  are  discrete  sums,  it 
is  clear  the  the  right-hand  side  of  (3.4)  is  a  r.v.  and  that  A^  e  F.  Observe 
now  that  if  X^ft)  is  right  continuous  with  left  limits,  then  the  same 
property  holds  for 


cV>  -  (  vt”  vt'  ’ 

^  0  ;  else 


Thus,  by  splitting  Xk  into  its  positive  and  negative  parts,  and  using  (3.4) , 
we  see  that  R^  is  a  r.v. 


We  will  use  Proposition  3.3  to  prove  that  R(n)  =  \  R,/n  is  a 

k=1 


randomized  I(n)  -  estimator.  Let  G=  o(X^  :  j  >  1)  (the  o-field  generated 


by  the  X j’s) • 

(3.5) .  Theorem.  Assume  A1  -  A6  and  A8.  Then  £{1*^  I  G}  =  1^. 

Proof .  Let  f  :  RP  R  be  a  non-negative  bounded  continuous  function  with 
compact  support.  Fix  an  arbitrary  selection  of  time  indices  t1f...,tp,  a 
write  f(Xk)  for  f (Xj^ti ) , . . . /xk(tp) ).  Assume,  for  the  moment  that  X^ 
non-negative .  Then , 

(3.6)  E{f(x"  )  r  x"  ( t )H  ( dt )} 

km  *  0  km  k 


CD 

-  I 

j=0 

CD 

-  I 

j-0 


E{f(x"  )X"  ( j2'n)(H.  ( j2~n)  -  H  ( ( j— 1 )2~n) ) } 

Km  kui  Jc  K 

E{f(x"  )x"  ( j2_n)(G( j2_n)  -  G(( j-1)2*n))} 

Km  Km 


=  E{f(x"  )  J"  X"  ( t)G(dt)} 
km  1 0  km 

where  the  second  equality  follows  from  independence  of  Xk  and  H^.  Since 
X^  is  right  continuous,  *  ^Xicm*  as  n  *  m  by  continuity  of  f, 

so  that  Proposition  3.3  implies  that 


(3.7)  li»  f(x"m)  i;  X^,tm  (it)  -  t(X  >  J”  V(t)Vdt)  .... 

n>* 


Now,  x£^  is  bounded  by  m  and  H^®)  is  integrable  (Lemma  3.1)  so 
evidently  the  left-hand  side  of  (3.7)  is  dominated  by  an  integrable  r.v. 
Hence,  one  can  apply  dominated  convergence  to  (3.6),  yielding 

B{f(X.  )  j"  X.  (t)H.  (dt)> 
km  1 0  km  k 


(3.8) 


=  E{f(X.  )  j"  X  (t)  G(dt) }  . 

km  J  0  km 


Since  f  has  compact  support,  f^X^)  =  ^xkM^  f°r  m  greater  than 

some  M.  Then,  the  non-negativity  of  f  implies  that  the  r.v.'s  on  both 


■ides  of  (3.8)  are  increasing  in  m,  for  m  >  M.  Hence,  monotone  convergence 
can  be  applied  to  (3.8),  yielding 

E{f(Xk)  J*  Xk(t)Hk(dt)} 

-  E{f(Xk)  J*  Xk(t)  G  ( dt )  }  . 

By  appealing  to  the  montone  class  theorem  ([7],  p.14)  and  using  the  fact  that 
the  functions  f(Xk)  generate  o(Xk),  one  obtains 

E{Z  j“  Xk(t)Hk(dt)}  -  E{Z  j“  Xk(t)G(dt)} 

for  all  bounded  Z  e  a(X^).  Hence,  by  definition  of  conditional  expectation, 

E{j;  Xk(t)Hk(dt)  I  Xk)  -  Jo  Xk(t)G(dt)  . 

In  particular,  E(J q |Xk( t ) | Hk<dt )}  ■  J  ^  E|Xk(t) |G(dt)  <  •  (see  A6),  so 
P(Ak)  *  0  for  k  >  1.  Thus, 

Jq  ^(tX^tdt)  -  \  a. s. 

,Hk>,  R^  is  independent  of 

o(Xj  :  j  /  k)  and  hence  ([6],  p.  308) 

|  G)  -  EtJ^  |  X^  «  IR  a.s.  , 

completing  the  proof  of  the  theorem  in  the  case  that  Xk  >  0.  For  the  general 
situation,  split  Xk  into  its  positive  and  negative  parts,  repeat  the  above 
argument,  and  recombine  using  the  fact  that  P(Ak)  =0.  | I 

(3.9).  Corollary.  Under  Al  -  A6  and  A8,  ERk  *  a. 

We  wish  to  prove  limit  theorems  in  terms  of  the  parameter  t,  where  t 
corresponds  to  an  index  of  effort.  In  the  applications  we  will  be  considering 
in  Sections  4  through  6,  the  effort  required  to  generate  Hk  will  be 


so  that  I  xk)  *  lk  a. 8.  Since  Rj^  e  o(Xk 


negligible  compared  to  the  effort  necessary  to  simulate  Ik.  Then,  sk(Tk) 
is  the  effort  required  to  generate  where  *  sup{t  :  H^t)  <  H^*)}. 

Note  that  and  are  F -measurable,  due  to  right  continuity  of 

and 


(3.10) .  Proposition.  Assume  A1  -  A6  and  A8.  Then 

(3.11)  E{Sk(Tk)>  Tk  <  •)  -  Jq  ESk(t)F(dt) 

where  F(t)  =  P{xk  <  tltEfS^T^;  A}  A  ES^X^I^,  where  IA(»)  is  1  or 
0  depending  on  whether  or  not  w  e  A). 

Proof.  Let  S^tt)  =  rain{Sk(t),m}  and  set  x£  -  2“n( j+1)  on 
{j  <  2ntk  <  j+1}.  Then, 

(3.12)  riS^lTj),  Tk  <  -) 


-  I  ECSk-|<T”),  2~"j  <  Tk  <  2-n(  j+1 )) 

00 

-  I  ES.  (2_n( j+1))P{2"nj  <  T  <  2~n( j+1 ) } 

j»0  kni  * 

where  the  second  inequality  is  due  to  independence  of  and  Tk.  Applying 

bounded  convergence  and  then  monotone  convergence  to  (3.12)  yields  (3.11).  || 

Let  M(t)  ■  max{k  s  S^x^)  +...+  Sk(xk)  <  t}»  then  M(t)  is  the  number 
of  Rk's  generated  with  effort  t.  The  next  lemma  shows  that  M(t)  is  no 
smaller  than  the  number  of  Ik's  generated. 

(3.13).  Lemma.  Assume  Al  -  A6  and  A8.  Then  M(t)  >  N(t)  a.s.  for  t  >  0. 
Proof.  We  need  to  show  Sk(xk)  <  S^T)  a.s.;  this  inequality  is  trivial  if 
T  ■  By  monotone  convergence, 

<•»)  -  Hj^  (T) }  *  lim  EfH^T+n)  -  ^(T)}  =  lim  G(T+n)  -  G(T)  =  0  . 


n-H® 


n+« 


The  non-negativity  of  H^*)  “  ^(T)  implies  that  H^*)  =  H^T)  a.s. 
proving  that  xk  <  T.  I  I 
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Let  R(t )  4  R(M(t) );  R(t )  is  a  randomized  estimator  constructable  from 


t  units  of  effort. 

(3.14).  Theorem.  Assume  A1  -  A6  and  A8.  If  p{xk  <  »}  ■  1,  then 
R(t)  ♦  a  a.s.  as  t  ♦  «•. 

Proof.  Since  ■  1  a.s.  (Lemma  3.1),  p{Hk(“)  >  0}  >  0,  so  that 

P{tk  >  0}  >  0.  Hence,  PtS^tT^)  >  0}  >  0  by  Proposition  3.10  and  A3.  This, 

in  turn,  implies  that  M(t)  <  *  a.s.  The  fact  that  <  «®  a.s.  forces 

S^Tj')  to  be  finite  a.s.  (see  Al),  which  assures  that  M(t)  ♦  •  a.s.  as 

t  ♦  “.  The  theorem  then  follows  from  the  strong  law  for  Tf(n)  and  Corollary 

3.  II 


To  obtain  analogs  of  theorems  2.3  and  2.8,  we  need  to  analyze 
2 

8(8^  ~  a)  .  Our  expression  will  involve  the  function  K(s,t),  where  K(s,t) 
EHk(s)Hk(t). 

A9.  EH^(t)  <  -  for  t  c  R+  . 

(3.15).  Lemma.  Under  A8  -  A9,  K(s,t)  is  the  distribution  function  of  a 
o -finite  measure  on  R2 . 

Proof.  Note  that 

Hk(s)Hk(t)  4  H*(t) 

for  s  <  t.  Hence,  K(s,t)  <  K(t,t)  <  •  so  K(s,t)  is  real-valued.  Let  sn 

tn  decrease  to  s,  t  respectively.  Since  H(s]c)H(t]C)  <  H(s.|)H(t1),  the 

finiteness  of  Kls^t^)  allows  application  of  dominated  convergence  to  prove 

that  K(s_,t  )  ♦  K(s,t);  K  is  therefore  continuous  from  above.  Also,  for 
n  n 

®1  ^  ®2'  ^1  ^  ^2# 


H(s2,t2) 


-  H(s1,t2)  -  H(82,t1) 


H(s1  ,t1 ) 


>  0 


so  that 


K(s2,t2)  -  K(s1  ,t2 )  -  Kfs^t,)  +  K(s1  ,t  1 )  >  0  . 

The  lemma  then  follows  from  Theorem  12.5  of  Billingsley  (1979). 


A10.  j 0  Jo  ®l*k(8>xk(t)|  K(ds,dt)  <•  . 

(3.16) .  Theorem .  under  A1  -  A6  and  A8  -  A10, 

(3.17)  ERjJ  =  J"  j"  EXk(s)Xk(t)  K(ds,dt)  <  •  . 

Proof.  For  X^  non-negative,  observe  that  for  any  integer  T, 

(3.18)  E{jJ  xjJm(s)xJn(t)H)t(ds)Hk(dt)} 

2"t  2°T 

-  I  l  E{x"  (j2“n)x"  (£2"n)AH  (j,4,n)} 
j-0  1=0 

2°T  2nT 

-  I  l  ex”  ( j2"n)x"  U2"n)AK(j,l,n) 

j-o  i-o  *■ 

*  Jo  Jo 

where 

AHk(j,l,n)  -  (Hk(j2”n)  -  Hk(( j-1)2'n))(Hk(12"n)  -  Hfc( (1-1 )2~n) ) 

AK(j,l,n)  -  X(j2'n,12"n)  -  K(  ( j-1 )2‘n,12"n)  -  K( j2-n, ( 1-1 )2_n) 

+  K(( j-1)2‘n,(l-1)2"n)  . 

Letting  n,  m,  and  T  tend  to  infinity  in  (3.18)  (in  that  order),  bounded  and 
monotone  convergence  proves  (3.17)  for  X^  non-negative.  For  X^  of  mixed 
sign,  split  X^  into  its  positive  and  negative  parts  and  recombine 
using  A10.  || 

Let  o2  -  J"  J“  EXk(s)Xk(t)K(ds,dt)  -  a2  -  Ed^-a)2  and  set  s(t)  = 

ES)'(t).  The  proof  of  the  following  theorem  is  identical  to  those  of  Theorems  2.3 
and  2.8. 

(3.19) .  Theorem.  Assume  Al  -  A6  and  A8  -  A10.  If  <  •)  ■  then 

lira  a(t)  •  |F(t)  -  a|  -  (a^  } .  s(t)F(dt))^2  a.s. 


(3.20) 


where  o"  <  •.  If,  in  addition,  ES^(Tk)  <  «,  then 


(3.21) 


t^  (R(t)  -  a)  ~ >  O  (J*  s(t)F(dt))1/^N(0,1 


38  t  ♦  •. 

2 

On  the  basis  of  the  CLT  (3.21),  it  is  natural  to  interpret  (^ES^T^/t  as 

the  asymptotic  variance  of  the  estimator  R(t). 

In  addition  to  the  estimator  R(t),  one  can  construct  a  second  estimator 

n 

based  on  the  framework  described  thus  far.  Let  h(n)  *  £  H.  (»)/n  and  set 

c  -  -  -  k"1 

^  I  R(n)A(n);  h(n)  >  0 

R(n)  =  ( 

^  0  ;  h(n)  m  0 

*+*  — 

The  estimator  R(n)  is  merely  R(n)  normalised  by  the  random  total  mass  of  the 
first  n  H^'s;  we  therefore  refer  to  R(n)  as  a  normalized  randomized 
I(n)  -  estimator,  and  set  R(t)  =  R(M(t) ) • 

The  following  result  has  an  identical  proof  to  that  of  Theorem  3.14  (recall 
that  EHk (•)  *  1  by  Lemma  3.1  so  that  h(n)  1  a.s.). 

(3.22).  Theorem.  Assume  Al  -  A6  and  A8.  If  P{tk  <  «*)  *  1,  then 
R(t)  ♦  o  a.s.  as  t  ♦  •. 

An  analog  to  Theorem  3.19  is  also  available,  under  a  certain  moment 
condition . 


AH.  j;  j;  E|Xk(s)X)c(t)|K(ds,dt)  <•  , 

A 

where  X^s)  “  (« >  ”  <* 


(3.23).  Theorem.  Assume  Al  -  A6,  A8,  A9,  and  All.  If  P{S.  (t.  )  <  •}  =  1,  the 


(3.24) 


lira  a(t)  •  I R( t)  -  a|  *  (0  J.  s(t)F(dt) 


)^  a.s. 


where  o*  -  J*  J*  EX^ (s)X^ (t)K(ds,dt)  ‘  ErJ  <  •  (R^  4  R)t  -o^).  If,  in 


addition,  ESk(Tk)  <  •,  then 


(3.25) 


t^(R(t)  -  a)  =***>  0  (J*  s(t)F(dt))1/2N(0,1) 


■v.  %/'vV\  */v"* 

/.  ,  /  V- 


as  t  ♦  •. 


Proof .  Note  that  if  X^Jt)  satisfies  A1  -  A7,  then  X^t)  also  does.  Under 

A 

All,  one  can  therefore  apply  Theorem  3.19  to  the  process  to  obtain  the 

relation 

(3.26)  Tim  a(t)  •  |  \  J*  Xk<s)Hk{ds)|  -  (0R  j"  s(  t)F(dt)  )^  a.s. 

ka 1 


But  observe  that 

.  M(t)  ^  .  *  M(t)  „ 

SitT  I,  Jo  V*’Vds>  ■  mo  J,  1 0  <V>  -  «>Vd*> 

k“ 1  k“ 1 

__  A 

-  h(M(t) ) (R( t)  -  O)  . 

Since  h(M(t))  ♦  1  a.s.  as  t  ♦  «•,  (3.24)  follows  from  (3.26).  A  similar  proof 

is  valid  for  (3.25).  I  I 


We  wish  to  emphasize  that  A10  and  All  are  not  equivalent  moment 
hypotheses.  Consider  a  case  in  which  Xj((t)  “  1  for  0  <  t  <  1  and  vanishes 
elsewhere,  with  0  <  G(1)  <  1.  Then,  A10  is  always  satisfied  but  All  is  valid 


only  if  K(*»,»)  <  •.  On  the  other  hand,  suppose  that  X^(t)  =  1  for  t  >  0. 
Then,  All  always  holds  but  A10  is  valid  only  if  K(«*,«»)  <  •  . 

*2 

As  in  the  case  of  R(t),  relation  (3.25)  suggests  that  oR  ESk(xk)/t  may 

A 

be  interpreted  as  the  asymptotic  variance  of  the  estimator  R(t).  The  goal,  in 

Sections  4  through  6,  will  be  to  determine  randomization  schemes  that  make 
2  *2 

or  ESk ( Tk )  and/or  0R  ESk(Tk)  as  small  as  possible. 

2  *2 

Certain  bounds  for  0  and  o_  are  available.  The  lower  bound 

R  R 

2  2 

o  >  o  can  be  obtained  from  Proposition  1.9,  under  A1  -  A8.  Theorem  3.5 
R 

asserts  that 

E{j"  Xk(t)Hk(dt)  |  G)  *  lk  -  a 

*2  2 

from  which  it  follows,  by  Proposition  1.9,  that  0R  >  o  •  Note  that  the  lower 
bounds  are  attained  for  Hk  *  G  a.s.  Upper  bounds  are  also  available,  under 
certain  conditions. 


•>.*  N.V  ■  .•  „*  -  *  .  • .  • .  • . 

w  *  .  *  w  •  « *  •  w  '  • '  •*.•«".* » ■  .  • .  • 

«  *  m  V.  ■  ,  »  •.  . 


’  v'  \  *  -  *  *  *  4 '*  *  .  ■  4  * 


.V •-  •. 


.»  ‘  •  a*  *V*  *  * 
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3 

(3.27) .  Proposition.  Assume  Al  -  A6  and  A8.  Then,  if  EH^(»)  <  », 

(3.28)  a*  <  j"  Ex£(t)G(dt)  •  EhJ(»)  -  a2 

(3.29)  02  <  (j*  EX2( t )G(dt)  -  a2)  •  EH2(«)  . 

The  upper  bounds  are  attained. 

Proof.  By  Cauchy-Schwartz , 

*k  *  Vi  Vt>Vdt))2  <  (J“  X2(t)Hk(dt))Hk(-)  . 

Replicating  the  approximation  argument  of  Theorem  3.16  proves  that 

(3.30)  ER^  <  J*  EX^(t)K(dt,«)  . 

Let  f  be  a  bounded  continuous  function  with  compact  support.  Then,  the 
definition  of  K  implies  that 

«• 

l  f( j2"n){K( j2_n,-)  -  K({ j-1)2'n,«)} 

j-o 

m 

-  e{  l  f(j2~n){H  (j2”n)  -  H.  (( j-1)2"n)}H  (-)}  . 

j-0  K  K  K 

Letting  n  ♦  m,  we  obtain 

(3.31)  j"  f ( t)K(dt ," )  -  E{j"  mjH^dtJH^-)}  . 

The  monotone  class  theorem  extends  (3.31)  to  include  all  bounded  Borel 
measurable  f.  But,  by  Cauchy-Schwartz, 

E{JA  Hk(dt)  •  Hk(»)>  =  E{J  A  H  (dt))^  •  (JA  Hk(dt)  •  Hk  («») 2  >  ^  > 

<  E{JA  Hk(dt)>  •  E{Ja  «k(dt)  •  Hk(«)2} 

<  JA  G(dt)  •  EHk(-) 

so  that  we  may  conclude  that  K(dt,«)  <  G( dt )EHk (•) .  Applying  this  inequality  to 


16- 


(3.30)  yields  (3.28).  Inequality  (3.29)  is  obtained  by  letting  X^(t)  play  the 
role  of  X^t)  in  (3.28). 

The  upper  bounds  are  attained  by  taking  H^ft)  =  0  for  t  <  and  1  for 

t  >  MR,  where  {M^}  is  an  i.i.d.  sequence  of  r.v.'s  having  common 
distribution  G.  I | 

The  inequalities  (3.28)  and  (3.29)  may  be  easily  adapted  to  provide  bounds 
on  the  moments  entering  hypotheses  A10  and  All,  respectively. 

(3.32) .  Proposition.  Assume  A1  -  A5  and  A8.  Then,  A6  -  A7  and  A9  -  All  are 

satisfied  if  EH^(-)  <  -  and  j"  Ex£(t)G(dt)  <  -. 

2 

Proof.  First,  observe  that  ElX^CtH  <  EX^(t)  so  A6  is  satisfied.  Also, 

EH^(»)  <  •  implies  that  X(“,“)  «■  EH^(“)  <  ",  yielding  A9.  As  for  A10,  the 
argument  of  Proposition  3.27  is  easily  adapted  to  show  that 

(3.33)  J*  E|Xk(s)Xk(t)|K(ds,dt)  <  J  J  EX2(  t)G(dt )  •  EHj^  (•)  <•  . 

A 

For  All,  substitute  X  for  X  in  (3.33),  and  observe  that 

j"  EX2(t)G(dt)  -  J*  EX2(t)G(dt)  -  o2  . 

2  2 
As  for  a7,  a  is  a  lower  bound  on  oR.  | | 

We  turn  now  to  a  variance  reduction  technique  that  is  sometimes  applicable 

to  the  randomized  estimators  developed  above.  Let  0  * 

tg  <  t^  <...<  tj^_1  <  t^  *  •  and  set 

«!«•>  -  P(\  '  '  1  Vi  <  'k  <  *i> 

for  i  ”  1,...,t.  Suppose  that  one  can  generate  independent  deviates  from  the 
measures  Q^,...,Q^  (see  Sections  5  and  6  for  examples)  and  that  p^  4 
p{ti_1  <  <  t^}  is  known,  for  1  <  i  <  t.  These  assumptions  will  allow  us  to 

stratify  our  sampling  scheme . 

Let  {m^  :  k  >  1}  be  a  sequence  of  integers  taking  values  in 
{l,...,i};  the  ra^ 1 s  define  a  sampling  order  for  the  Q^'s.  If  ra^  *  1,  one 
generates  from  (independently  of  the  (X^S^)  's) ,  and  sets 


•  k  e 


K  -  i’o  vtwi!,dt>  ■ 

If  *  {j  <  n  «  ■  1}  (1  <  i  <  £),  and  kn^  is  the  cardinality  of 


ttni#  one  sets 


S*(n)  1  I  P.  I  R1/kni 
i-1  1  jew  ,  3  ni 


\n)  -  l  p.  I  H®(- 
i-1  jew  .  3 


)  A  .  . 

ni 


Let  tJ  -  {t  :  H®(t)  <  H® (— >> ,  and  set  Ns(t)  - 

max{n  :  S^(T®)  +...+  Sn(t8)  <  t} ;  Ns(t)  is  the  number  of  R^'s  generated  in 

t  units  of  effort.  We  wish  to  study  the  stratified  estimators 
R*(t)  -RS(NS(t))  and  RS(t)  -  R®(N8(t)  )/h*l(N8(t) ) . 

A12.  knj/n  ♦  cA  >  0  as  n  ♦  •,  for  1  <  i  <  £. 

1.1  dt  -  E(Sk(Tk)  |  V,  <  tk<  tj}.  c\  -  e{l£  |  V,  <  ,k  <  tt)  - 

«<%  1  Vi  «  \  <  V2''  "i2  ’  I  Vi  «  %  <  V  - 

(EfR^  |  t^_j  <  <  t^}  ).  The  following  theorem  is  a  stratified  analog  of  the 

CLT's  (3.21)  and  (3.25). 

(3.32).  Theorem .  Assume  A1  -  A6  and  A8  -  A12.  If  ES]c(t]()  <  •*,  then 


(3.33) 


(3.34) 


t^2  (R8(t)  -  a)  — >  o^(  \  c.  d.  Y2  N(0, 1 )  . 

i-1 

t1/2(R8(t)  -  o)  — >  aRS(  l  c1di)1/2N(0,1)  . 


as  t  ♦  where  <4'  °L  <-  and  <4  -  l  p\o\/ V  4 


V  2  "2. 

1  PL  o  /c  . 

i-1  1  1  * 


"•  A  A  A  A  A  .*•  A 


•  A  ■  %  • 


clearly 


Proof.  Assumptions  A10  and  All,  together  with  <  », 

2*2  ^ 

guarantee  the  finiteness  of  0  ,  0  _,  and  \  c d  .  Under  A12, 

RS  R5  4  1  1 

1*1 


(3.35)  I  Sk(Tj)/n  >  I  t  I  S  <tJ)A  )  •  (k  /„)  .  I  e  d  .... 

k*1  i*1  jew  .  J  J  i*1 

ni 


as  n  ♦  «•.  But 


N  (t)  N  (t)+1 

(3.36)  l  8  (tr)/»  (t)  <  t/H  (t)  <  l  S  (T  )/»  (t> 

k=1  K  K  k-1  K  K 


on  {n® (t )  <  However,  as  in  the  proof  of  Theorem  2.2,  it  is  easily  argued 

that  Ns(t)  <  •  a. s.  and  Ns(t)  +  «•  a.s.  as  t  ♦  •.  Using  this  fact  in 

(3.36),  (3.35)  implies  that  t/Ns(t)  is  "squeezed"  between  two  terms  converg- 
l 

ing  to  \  c. d. .  Hence, 
i=1 


(3.37) 


NS(t)/t  ♦  (  \  c*  d, )  a.s. 
i*1 


The  central  limit  theorem  and  the  converging-together  lemma  ([2],  p.  25)  imply 
that 

(3.38)  n2  (  l  R®  -  ER®)A  ,  ~>  (O./c^  )N(0, 1 ) 


jew 


ni 


as  n  ♦  ••  Since  the  r.v.'s  \  R®  are  independent, 

J£“ni 


(3.39) 


ti2  (R® (n)  -  a)  -*>  (  \  P^o2./c  y*  N(0, 1 ) 

i*1  11  1 


as  n  ♦  Arguing  as  in  Theorem  2.8,  we  see  that  (3.39)  and  (3.37)  together 
yield  (3.33);  a  similar  proof  works  for  (3.34).  || 

We  now  turn  to  the  question  of  determining  the  c^'s. 

A13.  p^o*  >0  ,  1  <  i  <  t  . 

2  t 

(3.40).  Proposition.  Under  A13,  the  minimum  of  g(c)  A  0Rg(  \  c^d^) 


over 


i V% wi ivuth. v-atoh-vt; r. ■t-.t^v.  r; .».  rrr^ryrT7T7  VJ-  V-  “.'• 


>■ 

r&w 


$5 

ft 


ft 


'K 
«v 
>  ' 


Bs 

’•:< 


wu 

wo 

wo 


I 


a 


% 


A 


v*C* 


w> 


K 


{c  j  ^  >  0}  is  (  £  pi°idj^  )2 


and  is  achieved  at 


i=  1 


1/, 

c,  =  a  p.O./d/ 
~i  ~  i  i  i 


(3.41) 


o  =  I  P^./d'2 
i-1  11  1 


♦  A 2  /  (•  \  .♦ 
Similarly,  the  minimum  of  h(c)  A  0  (_  L  c  d  J  over  {c 

3  R5  #  __  M  1  1 7 

^  *  2  1a  2 

(  £  pi°idi  )  and  is  achieved  at 


i=  1 


>  0}  is 


i-1 


~  ~  1a 

ci  *  opioi/dr 


(3.42) 


o*I  piffi/di2 
i-1  1  1 


The  proof  involves  a  simple  application  of  Lagrange  multipliers;  see  Theorem 
3  of  Glynn  (1963)  for  a  similar  argument.  We  shall  now  show  that  stratification 
always  provides  a  variance  reduction. 

(3.43).  Proposition.  Assume  Al  -  A6,  A8  -  A12,  and  ES^x^)  <  *.  Then, 


(3.44) 


(3.45) 


(  i  Vi*f)2  ‘  “p'V’v 


i-1 


Proof.  By  Cauchy-Schwartz 


(  l  p ‘  (  l  p4»JK  j  Pi«i) 


i-1 


i-1 


i-1 


(3.46) 


=  H  Pi°JESk(V  • 

1=1 
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V  v  '  •**  .*• 

*  V*  4*  v  *  '  «  >  '  -  .  *  ,  . 


>  «■•  «-* 


Using  Cauchy-Schwartz  again  yields 


1  2  _ 2  1 


(3.47)  I  *  I  Pl«*^>  V,  «  \  ‘ 


i-1 


i-1 


>  “Sc  -  ( j,  V<V  Vi  <  Tk  <  VV 


i 


proving  (3.44);  (3.45)  is  proved  similarly.  || 

Utilization  of  the  constants  c^,  specified  by  (3.41)  and  (3.42) 

requires  a  "trial  run",  in  order  to  obtain  estimates  for  o^,  o^,  and  <v  If 
one  wishes  to  dispense  with  a  "trial  run"  and  is  willing  to  accept  a  sub-optimal 
choice  of  the  c^'s,  consider  using  ci  *  Pi*  Observe  that  if  c^  *  p^,  then 


t 

(  l 


i-i 


<  of  •  ES.  (T.  ) 
R  k  k 


; 


the  inequality  derives  from  (3.46)  and  (3.47).  Hence,  the  suboptimal  choice 

-  p£  yields  a  variance  reduction  over  a  nonstratified  sampling  plan. 

One  final  remark  is  in  order.  Suppose  that  there  exists  T  such  that  G(t) 

-  0  for  t  <  T  and  1  for  t  >  T.  Then,  an  argument  similar  to  that  used  in 
Lemma  3.13  proves  that  -  G  a.s. ,  in  which  case  R^  -  1^  a.s.  Hence,  our 
randomized  framework  yields  no  new  estimators  when  I  -  EX(T)  for  some  non¬ 


random  T 


Poiason  Process  Randomization 


In  this  section,  we  consider  a  class  of  randomized  time  integral  estimators 
for  which  the  H^'s  previously  defined  are  non-homogeneous  Poisson  processes. 
Specifically,  we  shall  assume  that 

(4.1)  Hk(t)  "  I 

for  some  X  >  0,  where  {n^  :  k  >  1}  is  a  sequence  of  i.i.d.  Poisson  processes 
for  which  EN^ft)  =  t>  of  course,  the  family  {N^  :  k  >  1}  is  assumed  to  be 
independent  of  ^Xk,Sk*  8  k  >  ^  satisfying  A1  -  A5.  The  following  proposition 
is  an  immediate  consequence  of  standard  properties  of  the  Poisson  process. 

(4.2) .  Proposition.  The  sequence  {H^  s  k  >  l),  as  defined  through  (4.1), 
satisfies  A8  and  A9.  Furthermore,  Eh£(—)  <  •  for  n  >  1. 

Calculation  of  the  distribution  K  is  easy. 

(4.3) .  Proposition.  For  0  <  s  <  t , 

K(s,t)  -  ^  G(s)  +  G(s)G( t)  . 

Proof.  The  independent  increments  property  of  the  Poisson  process  implies  that 
for  s  <  t, 

K(s,t)  =  E(N(XG(s))N(XG(t)) 

X2 

=  Ar  E(N(XG(s))(N(XG(s))  +  N(XG(t))  -N(XG(s))) 

X2 

=  — r-  EN2(XG(s))  +  (XG( s) (XG( t)  -  XG( s)  ) 

X2  x2 

=  |  G(s)  +  G(s)G(t) .  |  | 

Theorems  3.19  and  3.23  require  moment  hypotheses  on  the  process  X^. 

B 1 .  j"  EX2 (t)G(dt)  <  »  . 

Propositions  3.32  and  4.2  prove  that  Bl  is  sufficient  to  guarantee  A6  -  A7  and 
A10  -  All. 

(4.4) .  Theorem.  Under  Al  -  A5  and  Bl, 


T7T 


(4.5) 


(4.6) 


Proof .  Note  that 


ffR  =  °2  +  X  Jo  EXk(t>G<dt> 


°b  *  °2  +  X  <Jq  EX^t>G(dt>  '  °2) 


°2  =  E{J*  J*  Xk(s)Xk(t)G(ds)G(dt)>  -  a2 


j;  j*  EVs)Vt)G(ds)G(dt)  -  a2 

Jq  Jq  EX^siX^tJGfds^fdt)  . 


Hence,  the  result  follows  from  Proposition  4.3  and  the  fact  that  the 
distribution  Kj(s,t)  =  G(min(s,t))  induces  a  measure  on  R2  which  has  support 
on  the  diagonal  s  =  t.  | | 

Theorem  4.4  suggests  that  the  normalized  estimator  R(t)  has  better  large 
sample  properties  than  R(t)}  consequently,  we  shall  consider  only 

A 

R(t)  for  the  remainder  of  this  section.  We  shall  now  calculate  P. 

(4.7).  Proposition.  The  distribution  F(x)  =  P{tk  <  x}  is  given  by 


F( x )  =  exp(-X(1  -  G( x) ) ) 

Proof.  Evidently, 

(4.8)  Tk  =  sup{t  :  H^t)  <  H^"*)} 

=  inf {t  s  H^t)  >  Hk(-)} 

=  inf {t  :  Nj^UGJt))  >  Nk(X)} 

=  inf {t  j  G(t)  >  Tk(Nk(X))/X} 

where  Tk(j)  =  inf{t  :  ^(t)  >  j} .  Let  A(t)  =  P{t  -  Tk(Nk(t))  >  u} .  Then 

A( t)  solves  the  renewal  equation 


(4.9) 


where 


A(t)  =  a(t)  +  A(t-s)e  sds 


C  0  ,  t  <  u 

l  e  1  i  t  >  u 


Since  the  renewal  function  M{t)  =  EN^ft)  =  t  is  known,  (4.9)  can  be  solved 
explicitly  (see  [6],  p.  294): 

(  0  f  t  <  u 

(4.10)  A(t>  =  < 

^  e  t  t  >  u  . 

Now,  (4.8)  implies  that  {t^  <  x)  =  {Tk(Nk(X))/A  >  G(x)}  so  (4.10)  yields 

F(x)  =  P{A  -  Tk(Nk(X))  >  X(1  -  G( x) ) 

=  exp(-X(1  -  G(x) ) .  || 

We  shall  now  assume  that: 

B2.  s(t)  =  ESk ( t )  is  strictly  increasing. 

(4.11) .  Proposition.  Under  B2, 

ESk(Tk)  =  j“  {1  -  exp(-A(1  -  G(s_1 (x) ) ) )}dx  . 

The  proof  is  immediate  from  Propositions  3.10  and  4.5. 

Theorem  4.4  and  Proposition  4.11  can  be  used  to  optimize  the  value  of  X  so 
as  to  minimize  the  asymptotic  variance  of  R(t).  For  example,  suppose  that 
G(t)  =  t/T  for  t  <  T  and  1  for  t  >  T,  and  let  s(t)  =  t.  Then  the  goal  is 
to  find  X*  to  minimize  (see  Theorem  3.23) 

(o2  +  j  0)<jJ  {1  -  exp(-A ( 1  -  x/T ) ) } dx ) 

I  00  A  2 

where  P  =  JQ  EXk(t)G(dt);  differentiation  leads  to  a  non-linear  equation  which 
X*  must  solve. 

In  any  case,  the  inequality  exp(-Xu)  >  1  -  Xu  for  Xu  >  0  leads  to  the 

bound 

ES.  (T  )  <  A  f"  (1  -  G(s  1(x)))dx 
k  k  •'0 

Hence,  Poisson  randomization  with  intensity  X  is  more  efficient  than  I(n)  if 
(a2  +  B/X ) (X  J*  (1  -  G(s-1 (x)))dx)  <  02 s ( T )  . 
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5.  Point  Mass  Randomization 

Suppose  that  {(Xk,Sk)  s  k  >  1}  is  a  family  of  processes  satisfying  A1  - 
A5  and 

Cl.  J"  EX^(t)G(dt)  <  •  . 

Let  h  be  a  non-negative  function  for  which: 

C2.  {t  :  h(t)  >  0}  -  [y,»),  j“  h(t)G(dt)  -  1 

C3.  j"  Ij  (u)/h2(u)G(du)  <  •  . 

The  family  {h^  :  k  >  1}  of  point  mass  randomized  estimates  will  be  defined 
through  the  rule 

C  min{G(t),G(Y“)>»  t  <  ^ 

(5.1)  H.(t)  -  \ 

^  G(y-)  +  l/htM^j  t  > 

where  {M^  :  k  >  1}  is  a  sequence  of  i.i.d.  r.v.'s,  independent  of  the 
<Xk#Sk)‘s,  having  common  distribution  PtM^  e  dt}  *  h(t)G(dt).  With  the 
Hkas  defined  through  (5.1), 

\  -  Jq  I[0/f)(u)Xk(u)G(du)  +  X(Mk)/h(Mk)  . 

(5.2) .  Proposition.  Uhder  Cl  -  C3,  A6  -  All  hold. 

Proof.  First,  observe  that  Hk(t)  ■  G(t)  for  t  <  Y»  For  t  >  y.  Cl  proves 
that 

EH^t)  *  G(Y-)  +  B{l/h(\)l  <  t} 

*  G(Y")  +  J*  Ij  tj  (u)/h(u)  •  h(u)G(du) 

-  G(t) 

so  that  satisfies  A8.  For  A6  -  A7  and  A9  -  All,  we  apply  Proposition 

3.32.  It  is  evident  that  Hk(«)  *  G(y-)  +  1/Mm^)  so  that  it  suffices  to 
prove  that  EO/MM^))"*  <  <*.  But 
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Eh(M^)  -  J0  1^  ^(u)h(u)  J  •  h(u)G(du) 
which  is  finite  by  C3.  I  I 

It  should  be  clear  that  the  moment  assumptions  Cl  and  C3  are  considerably 
stronger  than  is  necessary  to  obtain  A10  and  All. 

(5.3) .  Theorem.  Under  Cl  -  C3, 

rG(s)G(t)  ;  min(s,t)  <  Y 

(5.4)  K(s,t)  =  <  G(Y-)(G(s)  +  G(t)  -  G(Y-))  +  J*  1  [Y,min(s,t)] (u)/h(u)G(du); 

^  min(s,t)  >  Y 

(5.5)  cl  -  e,txk>  »  62(*k)  -  a2 

(5.6)  ;2  -  sl(;k) .  82(x'k) 

(5.7)  EW  “  Jo  8(t)h(t)G(dt) 

where 

81(I)  -  lo  I[0,T)<»I[0,T)(t)K(.)X(t)G(d,)G(dt) 

+  2  J 0  Jo  I(0,T)",)IIt.-)(t>,,l*>,,(t)G(d,>Oldt) 

62(V)  *  J*  I(  (t)*X2(t>/h(t)G(dt)  . 

Proof.  Relation  (5.7)  is  immediate  from  the  definition  of  H^.  For  (5.4), 
observe  that  if  s  <  Y  <  t, 

EHk(s)Hk(t)  -  G(s)  •  EH^ (t )  =  G(s)G(t)  . 

If  Y  <  s  <  t,  then 

EHk(s)Hk(t) 

=  E{  (G(Y~)  +  XtM^/MM^I  [0  j  (Mk))Hk(t)} 

-  G(T->Mk(t)  *  E{^(Mk)/h2(Mk)I[0>8)(Hk)) 


+  E{X(Mk)A(Mk)I[0  j 


-  G(Y-)G(t)  +  j"  I[Y  j  (u)Ex£(u)/h(u)G(du) 

+  G(Y-)(G(8)  -  G(y->) 

2  “2 

which  is  (5.4).  The  expressions  for  0R  and  oR  are  obtained  from  Theorems 
3.16  and  3.23  by  appropriately  integrating  against  K(ds,dt)  as  defined  by 
(5.4),  and  by  using  the  fact  that  the  integrands  are  symmetric  in  (s,t).  | 

A  certain  amount  of  analytical  optimisation  can  be  performed  in  our 
current  setting.  Specifically,  one  can  determine  the  function  h  which 

A 

minimizes  the  asymptotic  variance  of  R(t)  over  all  point  mass  estimators  in 

A 

which  y  ■  0.  Note  that  if  Y  ■  0,  the  asymptotic  variance  is  v1(h)/t, 
where 

(5.8)  v1(h)  -  {J*  EXR(t)/h(t)G(dt)}(Jg  s(t)h( t >G(dt> >  . 

(5.9) .  Proposition.  v^(h)  >  (J Q(s(t)KXR( t) )^  G(dt) )2.  The  minimum  is 
attained  by 

(5.10)  h(t)  -  a(EXR(t>/s(t))1‘^  ,  a  >  0  . 

Proof.  Observe  that 

(j"  (smEXj^tn^Gldt))2 

-  (j"  (EX*(t)/h(t>  •  s(t)h(t))1/2G(dt))2 

<  J*  EX2(t)/h(t)G(dt)  •  J*  s(t)h(t )G(dt ) 

by  the  Cauchy-Schwartz  inequality.  To  check  that  h(t),  as  defined  by 

(5.10) ,  attains  the  minimum,  it  is  only  necessary  to  observe  that  h(t)  is 
undefined  only  on  a  set  of  zero  measure  with  respect  to  s(t)G(dt); 
similarly,  1/h(t)  is  undefined  nn  a  set  of  zero  measure  with  respect  to 
Ex£(t)G(dt).  || 
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* 

Hence ,  if  J0  (BXk(t)/*(t) )'2G(dt)  <  •,  the  density  h  defined  by 
(5.11)  h(t)  -  (EXk(t)/s(t ))^/j"  (ExJttJ/sUM^Gldt)  . 

A 

Minimizes  the  asymptotic  variance  of  R(t)  over  all  point  mass  estimators 
with  Y  ■  0.  Strictly  speaking,  (5.11)  may  not  be  in  the  class  of  randomized 
estimators  discussed  thus  far,  since  the  equality 

A 

E{ l/htM^ ) }  ^  <  t>  -  G(t) 

*2 

may  not  hold  due  to  the  possibility  of  EX^ft)  vanishing  on  the  support  of 
G.  If  that  occurs,  set  G(t)  ■  IA(t)G(dt )/jA  G(dt),  where  A  “ 

{t  :  EXj^(t)  >  0}.  Note  that  G  satisfies  A6  and  All  if  and  only  if  G  does; 
furthermore, 

A 

^  <  t}  -  G(t)  . 
a  -  J*  EXk(t)G(dt)  . 

*  ~ 

Thus ,  h  always  gives  rise  to  a  randomized  estimator  if  G  is  allowed  to 

play  the  role  of  G. 

A 

One  interesting  property  of  the  optimal  density  h  is  that  (5.11)  is 
independent  of  G,  modulo  the  normalization  constant.  Furthermore,  although 

A 

h(t)  can  be  estimated  via  a  trial  run,  that  may  be  undesirable.  In  such  a 

*2 

case,  it  may  be  reasonable  to  assume  that  EXk(t)  is  approximately  constant 
and  s(t)  linear,  leading  to  the  approximately  optimal  density 
t~  ^  /j"  t~  ^G(dt). 

A 

Analytical  information  about  the  optimal  density  h  is  available  also  in 
the  case  that  Y  is  positive.  The  idea  is  to  define  the  functionals 

J1 (h)  -  o*(h)(j"  th(t)G(dt) ) 

J2(h)  -  j"  h(t)G(dt)  , 

a2 

where  0R(h)  is  given  by  formula  (5.6).  We  wish  to  minimize  J^(h)  subject 
to  the  constraint  ^(h)  *  1  over  functions  h  having  support  [y,*). 
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Necessary  conditions  for  a  minimizing  h  may  be  obtained  using  variational 

* 

arguments  (see  Smith  [12]  ).  Formal  analysis  shows  that  a  minimizing  h 


should  be  of  the  form 


h(t)  -  o(EX2(t)/(s(t)  +  d))1^  , 


for  scow  constants  a  and  d,  when  y  >  0. 

Before  concluding  this  section,  we  offer  an  extended  version  of  the  point 
bums  estimator  studied  thus  far.  For  a  density  h  having  support  [y ,«•)  and 
satisfying  C3,  let  (H. .  ;  i  >  1,  k  >  1}  be  a  family  of  i.i.d.r.v.'s  for 


which  e  dt}  -  h(t)G(dt).  Set 


(5.12) 


Hik(t) 


min{G(t),  G(y-)}>  t  <  M, 


G(y-)  +  1/h(Mik);  t  >  Mik 


and  let  Hk(t)  *  J  H.k(t)/n.  In  the  presence  of  Cl,  A6  -  All  then  hold. 


.2,  .  12 


2  “2 


Let  <Jp  ( n ) ,  Op(n)  be  the  corresponding  values  of  oR,  o  for  the  current 


randomization  scheme. 


(5.13).  Proposition.  Under  the  assumptions  of  the  above  paragraph. 


(5.14) 


c2(n)  *  aj/n  +  (n-1)o2/n 


(5.15) 


o2(n)  -  c£/n  +  (n-1 )o2/n 


2  *  2 

where  o  ,  a  are  defined  by  (5.5)  and  (5.6). 

ft  K 

<»  n 
Proof.  Let  R^k  *  J0  Xj^fOH^tdt)  and  observe  that  R^  «  l  R^/n.  Then, 

the  exchangeability  of  {r,  .  :  i  >  1}  proves  that 


ERk  - +  ^  ERikRjk)/n2 


ERl. /n  +  (n-1  )ER„.  R_.  /n  . 


But  R 


has  the  sane  distribution  as  the  randomized  estimator  defined 


Ik 

2  2 

through  (S.  1);  hence  ER^k  ■  <JR*  To  calculate  the  expectation  of  the  product 
term,  note  that 


(5.16) 


R(RikR»  1  V  *  E(Rik'tRa  '  W  I  V 

■  RtRii.'(Ra  I  «k>  I  V 


1  V"R*  I  V  ’  xk  ' 


the  second  equality  follows  from  the  conditional  independence  of  0(Xk,H1k) 
and  cfX^H^)  given  a{Xfc)  ([6],  p.  308);  the  final  equality  was  actually 

proven  during  the  argument  of  theorem  3.5.  Relation  (5.16)  yields  ERi]<R2k  = 

2  * 
EI^,  which  proves  (5.14)»  (5.15)  is  obtained  from  (5.14)  by  substituting  X 

for  X.  |  | 


The  final  ingredient  in  calculating  the  asymptotic  variance  of  estimators 
based  on  (5.12)  is  the  value  of  ESR  ( Tk ) .  But  *  max{Mik  :  1  <  i  <  n}  ,  so 
that 


P{tk  <  x}  -  (J*  I[0  x](t)h(t)G(dt))n  t 

thus,  if  s(t)  is  strictly  increasing, 

Esk(Tk)  -  J-€j-  I  _1  ( t)h( t)G(dt) )ndx  . 

(s  (x),») 

It  should  be  clear  that  point  mass  randomization  is  amenable  to  variance 
reduction  via  stratification.  Observe  that  if  0  -  t0  <  T  <  t1  <  t2  <• • •< 

t  ■  •,  then  simulating  Mu  from  the  conditional  distribution 

"  fci 

h(u)Ir  .  . (u)G(du)/l.  h(u)G(du)  and  defining  Hu  via  the  rule  (5.1) 

lti_i,ti>  fci-1 

yields  a  deviate  from  Q^. 


a 


'.rV.V 


6.  Integrals  with  Random  Endpoints 


Assume  that  { (X. ,S. )  :  k  >  1}  Is  a  family  of  processes  satisfying  A1  - 


A5  and 


J*  EX^(t)G(dt)  <  • 


Let  L  be  a  probability  distribution  function  satisfying: 
D2.  L(O-)  =  0,  sup{t  :  L(t)  <  l)  -  T 


He  define 


(6.1) 


J*  G(du)/L(u-)  <  “,  where  L(t)  *  1  -  L(t) 
{H^  :  k  >  1}  by  the  rule 


Hk<t) 


j"  I  [0,t]  (u)G(du)  A.(u-)»  t  <  Nk 


J0  I[0fN  j  (u)G(du)/L(u-)»  t  >  Nk 


where  the  sequence  {n^  :  k  >  1}  is  a  family  of  i.i.d.r.v. 's,  independent  of 
the  (x]C*Sk),a#  such  that  Ptl^  e  dt>  -  L(dt).  With  the  l^’s  defined 


through  (6.1), 


*k~  Jo  I{0,llkJ(tt>Vtt>/E(tt')G<du) 


(6.2).  Proposition.  Under  Dl  -  D3,  A6  -  All  hold. 


Proof.  It  is  evident  that 


EHk(t)  "  J o  Jo  i[0,t](u)I(t,-)(8)A<u')G(du)L(d8) 
+  Jo  Jo  I[0,s](u)I[0,t](8)/i:(u')G(dU>I'(d8> 

“  j"  1  [0,t]  <u)L^t)  A*(u-)G(du) 

+  j"  1  [0,t]  "  I«(t) )  A(u-)G(du) 


*  G(t)  ; 

the  second  equality  is  Fubini's  theorem,  whereas  the  third  uses  the  fact  that 
L(t-)  >  0  on  [0,T] .  Hence,  satisfies  A8.  For  A6  -  A7  and  A9  -  All,  we 

apply  Proposition  3.32  -  it  will  suffice  to  prove  that  EH^f*)  <  ••  But 


•' 


SSfSRWi 


EH^(-)  -  E{j"  I[q  n  ](u)/L{u-)G(du)}3 

■  Ool'lo  ItO,*] 

“  ioiX  ^<®ax<u1»tt2,U3)“)  ^  1A<u  -)G(dui) 

i*1 

<  JqJq  S  1/L(u  -)G(du  ) 

i-1  1  1 

“  (J"  1/L(u-)G(du))2  <  -  .  || 


(6.3). 

Theorem. 

Assume  Dl  -  D3.  Then, 

(6.4) 

K(ds,dt)  *  L ( max (s,t )-)/£( s- )L(t- )G(ds)G(dt) 

(6.5) 

°R  "  r0  a(Xk;t)A(t-)G(dt)  -  a2 

(6.6) 

a2  =  j"  a(Xk>t)A(t-)G(dt) 

(6.7) 

ESk<Tk)  “  J"  s(t)L(dt> 

where 

a(Y>t)  - 

2  in  I/+  _ v ( S)EY( S )Y( t )G(ds)  +  EY2(t)(G(t)  -  G(t-))  . 

'0  (t,*) 

Proof. 

For  0  < 

s  <  t. 

K(s,t) 

“  eUo  i 

C0,Hk)(tt)I[0,s](o)^<,1')G(,S“>  '  Jo  I[0,Nk]<v>Il0.tl,v)^:,^)G<avl> 

- 

OoJo  I|0,x]l“>I[0,xl  <V)IC0,s]  ,o)I[0,t]  ('')A.(«-)L(v-)G(do)G(av)L(dx) 

- 

Mo  I[0, 

aj (u)I (v)L(max(u,v)-) /L(u- )L( v- )G(du)G(dv) 

proving 

(6.4). 

For  (6.5),  we  use  the  symmetry  of  EXk(s)Xk(t)  in  (s,t)  to 

write 

2 

ERk  * 

j"j"  EXk(s)Xk(t)L(max(s,t)-)/t(s-)L(t-)G(ds)G(dt) 

*  2  J  q  J*  (s)EXk(s)Xk(t)A(t-)G(ds)G(dt) 
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+  rX  I{t}(s)EX*(t)/L(t-)G(ds)G(dt) 

-  a(Xk»t)/L(t-)G(dt) 

A 

proving  (6.5);  for  (6.6),  apply  (6.5)  to  X^.  Since  !<(•)  ■  1,  <  •  a.s. 

so 

ESk(Tk>  ■  r0  s( t)F(dt)  -  j“  s(t)P{Tk  e  dt} 

-  J“  s(t)L(dt)  .  || 

As  in  Section  5,  one  can  try  to  analytically  optimize  the  choice  of  L. 
We  shall  require  several  additional  assumptions. 

D4.  s(t)  <  ■  for  t  >  0 

D5.  there  exists  a  non-negative  function  u(t)  such  that 

A 

a(Xk;t)/u(t)G(dt)  ”  s(dt)  (by  D4  and  A3,  s(t)  is  a  finite  non¬ 
decreasing  function) . 

A 

The  asymptotic  variance  associated  with  R(t),  for  a  fixed  L,  is  given  by 

A 

v2(L)/t  where 

VL)  "  X  »<«k»t)/E(t-)C(dt))(Jp  s(t)L(dt) )  . 

(6.8).  Proposition.  Under  D1  -  D5, 

v2(L)  >  (J*  atXk^J/u^  (t)G(dt))2  . 

Proof.  Observe  that 

(j“  a(Xkit)/u4  (t)G(dt))2 

-  (J“  (a^/tJAtt-))1^  (L(t-)a(Xk;t)/u^  (t)G(dt))2 
<  Jq  a(Xk;t)A(t-)G(dt)  •  j"  L(t-)a(Xk;t)/u(t)G(dt) 


by  Cauchy-Schwartz.  But 


JQ  L(t-)«(Xk»t)/u(t)G(dt)  -  J0JQ  I [t  -)(x)a(Xk»t)/u(t)L(dx)G{dt) 

-  J“J"  Ifn  . (t)a(X.  ;t)/u(t)G(dt)L{dx) 

•'O'O  (0,xj  k 


■  J"  s(x)L(dx) 


where  the  second  equality  is  by  Fubini,  and  the  final  equality  uses  D5. 

Assumption  05  requires  absolute  continuity  of  s(dt)  with  respect  to 

a 

G(dt)  on  [0, T]  and  non-negativity  of  a(Xk»t).  Hie  non-negativity  of 

A  * 

a(Xk ;t )  clearly  arises  when  Xk  has  non-negative  correlations  in  time. 
However  it  also  holds  under  other  conditions. 

A 

D6.  (Xk(t)  s  t  >  0}  is  a  stationary  process. 

(6.9).  Proposition.  Suppose  D1  and  D6  hold  for  G(t)  =  1  -  exp(-Xt), 

A 

X  >  0.  Then  a(Xk;t)  is  non-negative  and  decreasing  to  zero. 

Proof.  By  stationarity. 


(6.10) 


a(Xk»t) 


2  j"  EXk(0)Xk(s-t)Xe_Xsds 
2  X  J~  EXk(0)Xk(u)e”Xudu  •  e 


By  D1  and  Cauchy-Schwartz , 

E(X  Jp|Xk(s)|e"Xsds)2  <  j"  EXk(s)2Xe"Xsds  <  -  . 

Hence,  Fubini  proves  that 

0  <  B(j"  Xk(s)Xe'Xsds)2  "  EXk(s)Xk(t)X2e_X(s+t)dsdt 

-  2  j"j“  EXk(0)Xk(s-t)Xe"Xsds  Xe"Xtdt 

*>  X  •  j*  EX  (0)X  (u)e”Xudu  , 

*  0  k  k 


from  which  it  follows,  using  (6.10),  that  a(Xk?t)  is  non-negative  and 
decreasing  to  zero.  | | 


Assume  that  G(t)  *  e  Xt  and  that  s(t)  =  st.  Proposition  6.9  can  be 
used  to  find  a  distribution  L  which  is  approximately  optimal,  provided  that 

A 

is  reasonably  close  to  stationarity.  Substituting  the  expression  for 

A 

a(X^>t)  derived  in  the  proof  of  Proposition  6.9  into  the  defining  relation 
for  u(t),  one  gets 

(2X  f"  EX  (0 )X  ( t)e”Xudu)  •  e~2Xtdt  -  u(t)s  dt  t 
'  0  k  k 

_2\t 

hence  u(t)  *  8e  for  some  8.  The  following  property  holds  in  our 

current  setting. 

D7.  u(t)  is  a  continuous  function  which  decreases  to  zero  as  t  ♦  «*. 
(6.11).  Proposition.  Under  D1  -  D5  and  D7,  L^(t)  =  1  -  (u(t)/u(0) )^2  i8  a 

A 

distribution  function  which  minimizes  v^tL). 

Proof.  The  continuity  of  u  implies  that  L  has  the  required  right 
continuity,  so  that  L  is  a  distribution  function.  From  Proposition  6.8, 
V2(L)  >  U*  a<xk»t)/u^  (t)G(dt))2 

=  (J*  a(Xk»t)/L1(t)G(dt))(j"  a(Xk;t)L1(t)/u(t)G(dt)) 

=  (j“  a(Xk ;t )/tT1  (t)G(dt) )  ( JqJq  I[Q  (t)a(Xkft)/u(t)G(dt)L1(du)) 

=  (J*  a(Xk;t)/L1(t)G(dt))(j"  sCujL^du))  *  V^L^  , 

where  the  second  equality  is  by  Fubini,  and  the  third  follows  from  the 
defining  relation  for  u.  I  I 

As  a  consequence  of  Proposition  6.11,  we  see  that  L(dt)  =  Xe  Xtdt  is 

A 

approximately  optimal  for  a  simulation  problem  in  which  X  is  approximately 
stationary  and  G(dt)  =  Xe  Xtdt. 

We  now  turn  to  an  application  of  antithetic  variance  reduction  to 
integral  estimators  with  random  endpoints.  Let  L  ^ (x)  A  inf{u  s  L(u)  >  x} 
and  assume  that  there  exists  a  sequence  of  i.i.d.r.v.'s  {Uk  :  k  >  1}, 


independent  of  {(Xk,Sk)  :  k  >  1}.  Set  Nk1  =  L-1(Uk),  Nk2  =  L-1(1  -  Uk+1) 
for  k  odd;  Nki  =  L_1(1  -  Uk_.,),  N^  =  L_1(Uk)  for  k  even.  Let 

Nki 

Rki  =  JQ  Xk(s)/L(s-)G(ds) 

for  i  ”  1,2;  the  pairs  and  *Rk2,Rk+1  2*  odd)  are  said  to 

n 

be  antithetic.  Let  R«(n)  =  \  R  ./n.  The  following  theorem  involves  a 

k=1  * 

proof  similar  to  that  of  Theorem  3.19. 

(6.12) .  Theorem.  Assume  Dl  -  D4.  If  Jq  s(t)L(dt)  <  then 

(6.13)  t/2  (R^M^t))  -  a)  ==>  a2(Rn  +  R21)j;  s(t)L(dt)/2 

as  t  -*•  <*>,  where  M^(t)  =  max{2n  ;  +. ..+  N2n  1  <  t}. 

Based  on  Theorem  6.3,  the  antithetic  estimator  Rj (M^ (t ) )  is  more 
efficient  than  R(t)  if  o^(R^  +  R2^)/2  <  a2  =  02(R11). 

D8.  Xk(t)  >  0  -.s. 

2  2 

(6.14) .  Proposition .  Assume  Dl  -  D4  and  D8.  Then,  a  (R^  +  R21)  *  2  0R* 
Proof.  Clearly,  it  is  sufficient  to  prove  that  cov(R11,R21)  <0.  By  D8, 

L  1 (ui  )  _  L  1 (u2 )  _ 

R11  "  R22  =  J<)  Vs)/L(s-)G(ds)  -  JQ  Vs)/L(s-)G(ds) 

has  the  same  sign  as 


L*1(1-u2)  _  L_1(1-u1) 

r12  -  R21  =  Jq  X2(s)A(s-)G(ds)  -  }Q  X2(s)  A(s-*)G(ds)  . 

Hence, 

0  <  EfR^  -  R22)(R12  -  R21>  =  "2  cov(R11,R21)  .  II 

Variance  reduction  for  integrals  with  random  endpoints  can  also  be 

accomplished  via  stratification.  Note  that  if  0  =  t0  <  t,  <...<  tn  =  »,  the 

simulation  of  from  the  conditional  distribution  L(du)/(L(t^-)  -  L(ti_1~) 

s 

(tj_j  <  u  <  t^)  and  definition  of  ^  via  (6.1)  yields  a  deviate  from  Qi. 
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